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ABSTRACT 


In  this  part,  equations  for  pressure  distributions  and  generalized 
aerodynamic  forces  are  derived  for  a  thin  nonplanar  lifting  surface  in  simple 
harmonic  motion  at  subsonic  speeds.  A  digital  computer  program,  written  in 
Fortran  IV,  is  also  presented  herein.  The  computer  program  will  generate 
up  to  a  ten  by  ten  matrix  of  generalized  aerodynamic  forces  when  given  data 
for  the  geometry  of  a  planar  lifting  surface  with  a  folded  planar  tip,  the  flight 
Mach  number,  the  reduced  frequency  of  motion,  and  some  control  constants. 
Control  surface  deflections  are  not  accounted  for  in  this  study. 

The  kernel  function  method  given  by  Watkins,  Runyan,  and  Woolston 
(Reference  1),  which  relates  the  pressure  distribution  to  the  downwash  on  a 
planar  lifting  surface,  has  been  extended  and  applied  to  a  nonplanar  lifting 
surface.  Hsu's  technique  (Reference  4)  of  employing  Gaussian  quadrature 
formulas  is  used  when  integrating  the  product  of  the  kernal  function  and  the 
lift  function  over  the  planform  area. 

Recommendations  are  made  to  extend  the  method  to  account  for  blunted 
leading  edges  and  the  accompanying  airfoil  thickness  and  to  account  for 
control  surface  deflections. 
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I.  INTRODUCTION 


The  first  published  numerical  method  for  solving  the  subsonic  pressure 
distribution  problem  for  planar  lifting  surfaces  undergoing  simple  harmonic 
motion  was  developed  at  NASA's  Langley  Research  Center  by  Watkins, 
Runyan,  and  Woolston  (Reference  1).  Watkins,  et  aL  ,  presented  two 
methods  of  handling  the  numerical  integration  of  the  kernel  function  in  the 
region  where  high-order  singularities  exist.  Both  methods  involved  a  dense 
concentration  of  integration  points  in  the  neighborhood  of  the  singularity. 
Using  these  methods,  it  is  possible  to  obtain  downwash  integrals,  in  terms 
of  the  pres  sure -loading  coefficients,  at  any  arbitrary  set  of  points  on  the 
surface  (e.  g.  ,  at  all  the  kinematic  downwash  points  known  from  previously 
determined  vibration  mode  data).  However,  in  order  to  reduce  the  running 
time  on  the  computer,  the  downwash  integrals  were  obtained  at  a  selected 
set  of  collocation  points,  suchas  those  at  intersections  of  quarter,  half,  and 
three-quarter  chord  stations  and  like  half- span  stations.  When  downwashes 
were  matched  exactly  (and  thus,  boundary  conditions)  at  these  collocation 
points,  responsibility  was  placed  upon  the  user  to  evaluate  the  kinematic 
downwashes  there.  A  least-square  error  surface  fitted  to  the  mode  data  was 
commonly  used  to  evaluate  them.  Furthermore,  if  the  user  desired  that  the 
boundary  conditions  be  satisfied  at  a  greater  number  of  points,  it  was  neces¬ 
sary  that  he  use  a  correspondingly  greater  number  of  loading  functions. 

Procedures  were  then  described  by  Rodden  and  Revell  (Reference  2) 
and  the  correct  form  of  the  equations  were  presented  by  Fromrne  (Reference 
3)  for  calculating  pressure-loading  coefficients  which  match  a  greater 
number  of  kinematic  downwashes  than  coefficients,  in  the  sense  that  the  sum 
of  squares  of  amplitudes  of  differences  of  complex  numbers  are  minimized. 
Since  it  was  still  the  responsibility  of  the  user  to  evaluate  the  kinematic 
downwashes  at  the  collocation  point,  least- square  error  procedures  were 
used  twice:  once  implicitly  and  once  explicitly. 

Hsu  (Reference  4)  significantly  advanced  the  logical  development  of  the 
kernel  function  approach  when  he  established  an  optimum  set  of  collocation 
and  integration  points.  He  started  with  the  previously  established  chordwise 
pressure  functions  based  on  steady-state,  two-dimensional,  incompressible 
aerodynamics,  and  with  spanwise  loading  functions,  based  on  steady-state 
lifting-line  theory.  He  concluded  that  there  is  sufficient  reason  to  believe 
that  these  functions  display  the  proper  characteristics  near  the  edges  of 
lifting  surfaces  oscillating  in  a  compressible  fluid. 

Manuscript  released  by  authors  February  1965  for  publication  as  an  RTD  Technical  Documentary  Report. 
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Returning  to  the  two-dimensional  case,  Hsu  established  that  if  the 
chordwise  distribution  of  modal  deflections  (and  thus  downwashes)  is  accu¬ 
rately  represented  by  a  polynomial  of  degree  2N-1  and  is  approximated  by 
a  polynomial  of  degree  N-l,  then  the  integral  for  the  sectional  load  is  eval¬ 
uated  with  zero  error  by  a  N-point  Gaussian  quadrature  if  the  difference 
between  the  accurate  and  the  approximate  representation  of  the  downwashes 
is  made  equal  to  zero  at  each  of  the  N  points  (i.  e.  ,  the  chordwise  collocation 
stations).  Conversely,  still  for  the  two-dimensional  case,  Hsu  established 
that  if  the  product  of  the  pressure  function  and  the  kernel,  divided  by  the 
Jacobi-Gauss  weight  factor  (which  produces  the  square  root  singularity  at 
the  leading  edge),  is  accurately  represented  by  a  polynomial  of  degree 
2N-1,  then  the  integral  for  the  downwash  at  any  one  of  the  collocation 
stations  is  evaluated  with  zero  error  by  a  N-point  Gaussian  quadrature. 

These  N  points  are  then  made  the  chordwise  integration  stations. 

For  the  spanwise  direction,  using  lifting-line  theory,  Hsu  similarly 
established  M-spanwise  collocation  stations  and  M  +  1  interdigitated  spanwise 
integration  stations  plus  the  conditions  under  which  the  Gaussian  quadrature 
can  be  used  with  zero  error. 

It  is  important  to  note  that  the  kernel  of  the  integral  equation  for  the 
downwashes  in  unsteady,  three-dimensional,  compressible  flow  cannot  be 
accurately  represented  by  a  polynomial  of  finite  degree.  It  is  equally 
important  to  note,  however,  that,  because  of  the  edge  characteristics  of  the 
pressure  and  loading  functions,  the  Gaussian  quadratures  employed  at  Hsu's 
optimum  point  set  evaluate  the  integrals  with  the  least  squared  error  for 
a  given  number  of  integration  points.  We  have  yet  to  match  the  boundary 
conditions  using  Hsu's  method. 

The  downwash  matching  problem  in  Hsu's  approach  is  basically  the 
same  as  in  Watkin's  approach;  we  merely  have  a  more  logical  choice  of 
points  at  which  to  match  them.  In  the  examples  Hsu  used  to  demonstrate  his 
approach,  he  chose  to  use  the  same  number  of  pressure-loading  functions  as 
collocation  points.  However,  the  approach  is  not  dependent  upon  that  choice. 
If  a  smaller  number  of  pressure-loading  functions  are  used,  then  the  pro¬ 
cedures  described  by  Rodden,  Rcvell,  and  Fromme  may  be  used  to  compute 
pressure-loading  coefficients  which  yield  a  minimum  sum  of  squares  of 
amplitudes  of  differences  in  downwashes. 


A  need  has  arisen  for  application  of  the  kernel  function  method  to 
rion-planar  lifting  surfaces  on  future  aerospace  vehicles.  Application 
is  also  required  to  more  conventional  non-planar  surfaces  such  as  T-tail, 
V=-ta!L?L^  and  wing-vertical  tail  combinations. 

Professor  H.  Ashley  outlined  the  application  to  the  folded  tip  con¬ 
figuration.  A  computer  program  based  on  Ashley's  work  was  developed  for 
steady-state  flow  by  L.  Johnson,  et  al.,  of  the  Los  Angeles  Division  of 
North  American  Aviation,  Inc. 

The  work  reported  herein  is  based  on  Professor  Ashley's  outline. 
However,  the  expression  for  the  kernel  has  been  greatly  simplified  by 
Dr.  K.  R.  Rodemich  of  North  American  Aviation,  Inc.,  Space  and  Information 
Systems  Division. 
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II.  FUNDAMENTAL  EQUATIONS  OF  FLUID  MOTION 


Consider  a  body  immersed  in  a  compressible,  nonviscous,  perfect 
.fluid  and  assume  the  fluid  flow  to  be  isentropic  and  irrotational.  Under  these 
conditions,  a  velocity  potential  exists: 

q  =  (1) 

where  q  is  the  velocity  vector  of  a  fluid  element  and  V  is  the  gradient 
operator  (See  Reference  4.)  Also  under  these  conditions,  the  isentropic 
(constant  entropy)  pressure-density  relationship  is  valid.  Thus, 

a2  *  (!)*  = Y  RT  <2> 

Other  equations  which  govern  the  flow  are  the  continuity  equation  for 
conservation  of  mass 

!^+v(pq)  =  o  (3) 

and  Euler's  equations  for  conservation  of  momentum 


Dq  -1  n 
Dt  =  7  p 


where,  — ,  the  substantial  derivative,  is 


JD 

Dt 


(5) 


These  equations  may  be  combined,  as  described  in  Reference  5,  to  yield 
the  nonlinear,  unsteady  flow  equation 
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Consider,  then,  that  the  fluid  motion  consists  of  a  perturbation  super¬ 
imposed  on  a  uniform  stream  velocity  V=  Ui  +  Wk  parallel  to  the  xz-plane  •> 

of  a  rectangular  Cartesian  coordinate  system.  Then  the  velocity  potential 
may  be  expressed  as  the  sum  of  a  uniform  part  and  a  perturbation  part 

4  =  Ux  +  Wz  +  Cp  (7) 

and,  similarly,  the  velocity  vector  becomes 

q  =  V  +  =  V<p  (8) 


The  pressure  coefficient  at  any  point  in  an  isentropic  flow  field  is 


C 

P 


P  -  Pco 

1  ^T2 

2  p“  U® 


(9) 


where  p  -  p^  is  the  difference  between  local  pressure  and  free -stream  pres¬ 
sure,  Um  =  |V1  ,  and  1/2  n_  IU2  is  th«  Ay~Z~iz  Fiuin 

Kelvin's  equation  {Reference  5)  for  isentropic  flow 


Cp  =  2 


1  + 


(1  - 


_  _  ,  _  a  <p 
q*  q+2gr 


u: 


Y 

Y-l 


(10) 


A  complete  statement  of  the  fundamental  problem  requires  specifica¬ 
tion  of  the  boundary  conditions.  The  boundary  conditions  at  infinity  depend 
upon  the  free-stream  velocity.  When  it  is  less  than  the  speed  of  sound  in 
the  fluid,  the  disturbances  to  the  flow  die  out  and  are  not  felt  at  infinity. 

When  it  is  greater  than  the  sonic  speed,  then  in  the  region  where  disturbances 
are  felt,  even  at  infinity,  the  component  of  flow  due  to  the  disturbance  is 
directed  away  from  the  source  of  disturbance  and  otherwise  the  free-stream 
flow  is  undisturbed.  The  boundary  conditions  at  the  surface  of  the  body 
require  that  the  flow  be  tangent  to  the  surface  everywhere  on  the  body.  This 
condition  is  satisfied  by  the  equation 


D_ 

Dt 


£  (x,  y,  z,  t)  =  0 


(11) 
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where 


B(x,  y,  z,  t)  =  0  (12) 

is  the  equation  for  the  position  of  the  surface  at  any  time  t,  and  the  substan¬ 
tial  derivative  D/Dt  is  defined  by  Equation  5. 


III.  LINEARIZED  EQUATIONS  OF  MOTION 


Linearization  of  the  equations  of  motion  is  not  dependent  upon  an 
explicit  form  of  the  body  equation,  Equation  12,  so  long  as  the  normal 
derivatives  of  the  equation  are  everywhere  nearly  perpendicular  to  the  free- 
stream  direction.  Thin  lifting  surfaces  at  small  angle  of  attack  satisfy  this 
condition  and  are  treated  herein  and  in  Parts  2  and  4  of  this  report.  The 
special  considerations  required  for  thick  bodies  and  high  angles  of  attack  are 
treated  in  Parts  3  and  5.  The  following  development  is,  therefore,  restricted 
to  thin  airfoils. 

We  first  obtain  the  specialized  form  of  Equation  7  when  the  uniform 
stream  velocity  lies  along  the  x-axis;  i.  e.  ,  W  =  0  and,  therefore,  V  =  Ui, 

Uo,  =  U.  The  velocity  potential  is 

4>  =  Ux  +  <P  (13) 

and  the  velocity  vector  of  a  fluid  element  becomes 

q  =  (U  +  u)  i  +  vj  +  wk  (14) 


where 


u 


a  <p  b<p  ,  a  <p 

v=-,  andw=- 


The  perturbation  velocities  u,  v,  and  w  are  assumed  to  be  much  smaller 
than  the  free -stream  velocity;  i.  e. ,  u,  v,  w  «  U. 

The  linearization  procedure  when  applied  to  Equation  10  yields  the  fully 
linearized  pressure  coefficient 


00 


_0_ 

at 


+  Ur-)  ^ 

ax' 


(15) 
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and  when  applied  to  Equation  6  yields  the  fully  linearized  unsteady  flow 
equation 


8 Z<p  82  <P  8 Z<P 
8x2  8y2  8z2 


”■7 *•>»&*£ 


=  0 


(16) 


Next,  we  write  the  body  equation  for  a  thin,  nonplanar  lifting  surface 
(Figure  1),  in  terms  of  a  curvilinear  coordinate 


s  =  s  (y) 


s  represents  the  integral  of  distance  along  the  line  of  the  mean  position  of 
the  airfoil  from  the  centerline  to  y, 


s 


ds 


(17) 


The  position  of  the  surface  in  terms  of  s  and  n  {the  normal  to  S),  is 
separated  into  two  parts;  one  part  for  the  upper  (or  inner)  surface,  and  the 
other  for  the  lower  (or  outer)  surface 


B  (x,  n,  s,  t)  =  n  -  n  (x,  s)  -  n  (x,  s,  t) 

U  '  T 


(18a) 


B-  (x,  n,  s 


t)  = 


/v 

9 


/  v 


«  4-\ 

w9  V/ 


M 

\AOWi 


where  nT(x,  s)  represents  the  thickness  of  the  airfoil  and  n(x,  s,  t)  repre 
sents  the  elastic  deflection  of  the  airfoil.  In  accordance  with  Equation  11, 
we  use  the  operator 


_  •  9  .  •  '  ^  .  i*  9 

7  =  i  - —  +  l  —  +  k  - — 
8x  J  3s  8n 


(19) 


on  Equation  18a  to  get 


VB  (x,  n,  s,  t'/  =  -i 
u 


h  nT(x>  B)  +  h n(x-  ‘>1 

j 


8x 


.  / 
-J 


—  n.  (x,  s)  +  —  n(x,  8,  t)  +  k' 


Substitution  into 


8B 

mBu(Xl  s-  "•  *  IT  +  (ui  +  ’  7Bu=0 


(20) 


of  the  equation  for  the  upper  surface,  after  higher  order  terms  have  been 
discarded,  gives 


8  <? 
8  n 


)  n(x,  s. 


t)  +U  — n  (x, 
9x  t' 


b) 


(21) 
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The  same  procedure,  for  the  lower  surface,  gives 


9  P 
9n 


iL 

'9t 


+  IT- 


9x 


)  n  (x,  s;  t)  -  tj—  nT(x,  s) 


(22) 


Finally,  we  restrict  the  analysis  to  that  class  of  problems  in  which  the 
effects  of  thickness  on  the  time  dependent  forces  can  be  neglected.  By  letting 


nT(x,  s)  =  0 

we  get  a  single  expression  for  the  boundary  condition 


(23) 


9?  .9  -I  9  \  f  M 

9xT=  9t  *  U  9x  n  X*  *> 


(24) 


It  is  evident  from  Equation  24  that,  when  the  motion  of  the  surface  is 
simple  harmonic  motion, 


n(x,  8,  t)  =  n  (x,  s)  e 


iut 


(25) 


then, 


P  (x,  s,  n,  t)  =  V(x,  s,  n)  e 


iwt 


(26) 


Substitution  of  Equations  25  and  26  into  Equations  15,  16,  and  24 
give  8 


DP 

Dt 


2- 
V  P 


9  P  _  Dn 
9n  Dt 


(27) 


(28) 


(29) 
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where 


_d_ 

8x 


+  iw 


(30) 


and 


2 

V 


(31) 


To  this  degree  of  approximation,  M  =  and  a  =  affl. 


IV.  THE  ACCELERATION  POTENTIAL 


PLANAR  WINGS 


Equation  27  shows  that  the*,  calculation  of  pressure  requires  taking 
derivatives  of  the  velocity  potential,  and  Equation  29  states  the  boundary 
condition  which  must  be  satisfied.  Watkins,  Runyan  and  Woolston 
(Reference  1)  solved  this  problem  for  planar  surfaces  in  terms  of  a  series 
of  pressure  functions,  or  acceleration  potential  functions  (»};),  where, 


«1»  (x»yi  z)  =  ~  <f  (x,y,  z) 


(32) 


Integration  of  Equation  32  gives  the  general  expression  for  the  velocity 
potential  in  terms  of  the  acceleration  potential: 


(\,y,z)d\  (33) 


.  x 
-VA  ~ 


4>  (x,y,  z)  =  —  e 


5/ 

*  -a 


10) 


U  4> 


if  the  velocity  at  infinity  is  V  =  Ui.  The  velocity  potential  4>  due  to  a  pulsating 
doublet  satisfies  Equation  28  (in  the  planar  case  s  =  y  and  n  =  z);  and,  because 
the  order  of  operators  is  interchangeable,  the  acceleration  potential  ip  also 
satisfies  Equation  28. 

A  complete  discussion  of  the  application  of  boundary  conditions  in  the 
planar  case  is  given  in  Section  6-4  of  Reference  5.  The  application  in  the 
nonplanar  case  is  discussed  in  less  detail  in  the  following  text. 

NONPLANAR  WINGS 

In  the  nonplanar  case,  the  acceleration  potential  at  a  point  (x,  *,  N) 
due  to  a  pulsating  doublet  located  at  the  point  ((,  <■,  «),  or  (£ ,  q ,  l, )  in  the 
direction  of  n  is 


(x,  s,N) 


(34) 


where 


9n 


8  8 
cosy  (tj)  —  -  siny  (n) 


(35) 


r  =  V(x  -  6  )2  +  fl2  [(y  -  n)2  +  (z  -  £)2]  (36) 

and  Y(q)  is  angle  between  the  wing  and  the  xy-plane  at  the  point  { £,  q ,  £ ). 
The  perturbation  velocity  potential  may  be  built  up  from  a  distribution  of 
doublets  of  acceleration  potential  over  the  wing._  If  A  is  the  infinitesimal 
doublet  strength  at  (£,<r  ,q  ),  the  contribution  to <f  from  the  doublet  at  this 
point,  from  Equations  33  and  34,  is 


1  Ci) 


-  -A  9 

A#x,s,N)  =  17  9n  e 


X  -  g  X  -  £ 

-1GJ 


\  MX  R1 

U  +  a2  2 
a  0  a_B 


~  / 


GO 


-  CO 


where 


r*  =  V^2+  62  hy-^)2+  (z  -  &>2] 

and  the  velocity  component  normal  to  the  surface  at  (x,  s,N)  is 


’d\ 


(37) 


lim  9  — ,  VT, 

Aw{x,  s,  0)  =  N_0  —  Af(x,s,N) 


(38) 


where 


9 

9  N 


cosy  (y)  -r- 
o  z 


-  siny  (y) 


(39) 


Note  that  when  the  operator  9/9n  is  applied  in  Equation  37,  the  partials 
9/9 £  and  9/9q  maybe  replaced  by  -  9/9z  and  -  9/9y,  respectively. 
Substitution  of  Equation  37  into  Equation  38  gives 


Aw(x,  s,  0) 


A. 

U 


-  i  to 


e 


-x  zA 
u 


lim 
N  —  0 


/ 


X"  ^  eico(X  -  MR')/U02 


d\ 


(40) 


where  the  operator  P  is 
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P=  [cosY(y)  -  sinY(y)-^r]  [c08Y(n)~  -  »inY(n)~)  (41) 


and  finally 


A  _  Apd£dcr 


and  Ap  is  the  complex  amplitude  of  the  difference  in  pressure  on  the  upper 
and  lower  sides  of  the  surface  at  (£  ,c r), 

Ap (§  , cr )  =  pu  (£,<r)  -  P|  (£,cr)  (43) 

and  d£dcr  is  the  incremental  area  of  the  doublet  sheet. 

The  normal  wash  at  (x,  s, 0)  given  by  Equation  40  is  that  due  to  a 
point  pressure  doublet  at  (£  ,  <r  ,  n  =  0).  The  total  normal  wash  is  the  integral 
over  the  surface  of  all  the  pressure  doublets, 

ff  AP  (£  »or)  K(x  -  §  »  s,  <r ,  w,  M)d£dor  (44) 

J  4t r  pU 

where  j:  denotes  the  Mangier  formula  for  evaluating  infinite  integrals 

(Reference  7  ),  and  the  kernel  of  the  integral  equation  is  (omitting  the 
arguments  u,  M  for  brevity) 


ux 

o 

\  lim  1  u  p 

K(VSl,r)  =  n-0  e  P 

N  — ■  0 


0  iw(\-MR')/UB 
e 


d\  \  (45) 


where 


R«  =  WX2  +  B2r2 


/  2x  2 

r,  .=  -v  /  y  +  z 
1  V  o  o 


x  =  x  -  £  ,  y  =y-rl»  and  z  =  z  -  £ 
o  *  1  o  1  o 
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Now,  by  putting  k,  =  ajri/U  and  v  =  X/  Sr  ;  then,  by  putting 

/  2”  1 
v  -  M  vi  +  v 

u  = - - - the  integral  in  Equation  45  may  be  written 


r  °  ico(X  -  MR')/U62  r®  "lklU 

io=-/  - - 5T - dX  =  /  7=f- 

-CD  U^  VI  +  u 


where 


x  -  MR 
o 

Ul  = 

1 

By  breaking  up  the  interval  of  integration  into  three  subintervals  and 
in  the  first  two  integrals  letting  u  =  w/i 


1  -k,w 


r  l  r  1 

=  -1  /  6/  -•  -3-  dw  +  /  Vx 


0  V  i  «  w 


1  v  w  -  1 


1  -iku 
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0  V  i  +  u 


1  -k.w 


I  =  K  (k 
o  o 


l1  '  1  /  X  "7  dW  ‘  / 


1  -ik^u 
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du  (47) 


0  Vi  .  w 


0  V  i  +  u 


where  K0(kj)  is  the  modified  Bessel  function  of  the  second  kind  and  zeroth 
order  of  the  argument  kj. 

To  obtain  the  analytic  form  of  the  kernel,  we  write  the  operator  P 
(Equation  41)  in  a  more  convenient  form,  taking  advantage  of  the  fact  that 
I0  is  a  function  of  only  r^  when  xQ  is  held  constant 


PI  = 
o 


=  cos  [v(y)  -  Y(n>]  ^7“ 


+  [zoCOSY(y)-yo8inv(y)]  [Vosy(n)  -  y^invU,)]  (*“  )  (  ~  "Fr^" ) 
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The  resulting  kernel  is 


where 


K(xo,s,<r) 


2 

STIP 


K{xqIs  ,  tr ) 


X' 


TlKl(xo»s><r)+  T2K2(xo»8»(r) 
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(48) 


T,  =  cos  [  v(s)  -  v(ff)] 


To  = 


Z  y 

-O  .  .  v  Zo 

~  cosv(s)  -  —  sinv(s) 


Z  y 

-o  Lo 

—  cosy(£)  8inY{£) 


(48a) 
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dw 
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du  +  k 


2  f  iL_S _  , 

j  J  ) - -  du 


/ - ^  UU  T  K.  I  — 
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(48b) 
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U1  =  - 


X  -  MR 
o 


02r. 


(48c) 


and  K^(kj)  and  are  modified  Bessel  functions  of  the  second  kind  and 

first  and  second  orders.  The  sub-bar  indicates  division  by  Srj.jp,  e*g. , 
rl  =  ri/sTIP. 


EO 


V.  THE  BOUNDARY  CONDITIONS 


The  remainder  of  the  problem  is  tc  match  the  boundary  conditions;  i.  e., 
to  find  a  pressure-loading  function  Ap(£,  c r)  which,  when  inserted  into 
Equation  45  and  integrated  over  the  surface,  yields  the  kinematic  down- 
washes  at  selected  points  on  the  surface  w(x,  s,  0). 

In  subsonic  flow,  the  behavior  of  the  pressure  distribution  is  known  in 
the  area  of  the  wing  edges  from  a  few  of  the  exact  solutions  in  lifting  surface 
theory.  In  the  neighborhood  of  the  leading  edge,  the  pressure  should  behave 
as 

lim 


6—0 

In  the  neighborhood  of  the  trailing  edge  and  all  edges  parallel  to  the 
free-stream  direction,  the  pressure  should  behave  as 

lim  y/T 

6—0 

where  6  is  the  distance  to  the  wing  edge.  Both  Hsu  and  Watkins  employ  a 
linear  superposition  of  functions  that  satisfy  these  conditions.  Hsu's 
function  differs  from  Watkins  only  in  that  for  any  given  number  of  terms  in 
the  series,  Hsu's  terms  are  linear  combinations  of  Watkins  terms.  We  use 
a  normalized  form  of  the  function  given  by  Watkins: 


a  <rmf  (£ )  (49) 

nm  -  n  '  ' 

where  _ 

yi ) ■ ^/^f 

i  +  i 

f  (f)=  \/l  -I2U  (£);  1  <n 


Ap{£  ,cr)  =  -■ 


,u2/2 


N 


b(<r) 


M 


I 

n=0  m=0 


The  anm's  are  unknown  pressure  coefficients  to  be  determined  by  matching 
the  kinematic  downwashes  at  the  selected  points  (xj,  sr)  on  the  surface. 
Substitution  of  Equation  49  into  Equation  44  leads  to  the  matrix  equation 
given  by  Rodden  and  Revell  (Reference  2)  Equation  39,  for  the  point  set 
Xj,  sr 


We  now  reexamine  the  fundamentals  of  the  problem  before  proceeding 
to  evaluate  the  integrals  in  Equation  51  and  thence  to  solve  Equation  50. 

One  of  the  basic  reasons  for  development  of  the  kernel  function  method 
is  that  pressure  distributions  over  a  continuous  lifting  surface  are  smooth 
continuous  functions  that  can  be  represented  with  reasonable  eccuracy  by  a 
series  of  analytic  functions.  We  point  out  that  fn(£  )  can  be  written 


u  (I);  1  >  n 
n 
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and,  therefore,  the  inner  integral  in  Equation  51  may  be  written 


!1  ‘ 


■  / 


1  ) 

— - d£ 


(52) 


-1  >/l  -  I2 


where 


P  (l)  =  (1  -I)K(x.  -  S.s  .£) 

o  J 

p  (I)  (i  - 1 2)u  (l)K(x.  -  e,sr,£>;  i 

n  n  J 

Now,  we  assume  for  the  moment  that  the  kernel  K(2G  “  8r*<r^ 

represented  with  reasonable  accuracy  by  a  polynomial  in  £  .  Then,  Pn(£  ) 
is  also  a  polynomial  in  £,  and  the  Chebyshev- Gauss  quadrature  formula  may 
be  used  to  obtain  the  exact  value  of  the  integral  expression  (52);  1.  e. , 

1  P  (I)  -  V  ir 

(  ~^=^I  i  W+E 

J  k-‘ 

where 


(52a) 


E  = 


2tr 


22K(2K)! 


p  (2K)(X) 

n 


and, 

|X|  <1.0 

The  error  term  E  is  zero  if  Pn  is  a  polynomial  of  degree  <  2K  -  1. 

A  more  accurate  formula  which  utilizes  the  fact  that  Fn(1.0)  =  0  is 
used  by  Hsu  (and  by  us) 


Vs  1 


yi> 


(l  +£)  un(|) 
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Then,  the  inner  integral  in  Equation  51  may  be  written 


!1  =  /  V1- T_I Fn 

-1  v  1  +  £ 


U)d| 


(53) 


where 

F  (I)  =  K(x.  -  e,  8  ,  cr) 
o  J  -r  — 

F n(i)  =  (1  +1)  Un(f)  K  (x.  -  e.  8r  1  <n 

This  is  evaluated  by  the  L- point  Jacobi- Gauss  quadrature  with  the  weight 
function  yj  (1-  -  f ) / (1  +  %  )  (see  Reference  8,  Chapter  8).  The  resulting 
formula 


Is 


L 


2 


W.F  (SJ 
k  n  k 


(54) 


is  exact  if  Fn(£)  is  a  polynomial  of  degree  <  2L-1,  which  corresponds  to 
degree  2L  for  Pn(I ).  Putting  |  =  -cos  0,  the  polynomials 


are  orthogonal  with  respect  to  the  weight  function  V(1  -§)/(!  +  £ ) 


1 
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+jf  )♦„<!) 

m  n 


0,  m  £  n 


r» 


m  =  n 


This  is  easily  verified  by  expressing  the  integral  in  terms  of  0.  Referring 
to  formulas  in  Reference  8,  Chapter  8,  Section  8.4,  it  can  be  shown,  using 
these  polynomials,  that  Equation  54  takes  the  form 


N 

-1 


i  +  i 


=  F„«>d«  ■  lufi 


L 

2 
k  =  1 


H,  Uk)  (55) 
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where 


Hk  = 


and 

r  /2k  -  1  \ 

t.  =  -cos  — - r  irl 

*  k  \2L  +  1  / 

In  two-dimensional,  steady,  incompressible  flow,  there  is  an  optimum 
set  of  chordwise  collocation  stations  (xi)  for  the  determination  of  sectional 
lift,  depending  upon  the  order  of  the  polynomial  required  to  adequately  repre¬ 
sent  the  downwash  distribution 


x. 

J 


-cos 


r  2j 

2N  +  1 


1,  2, 


N. 


(56) 


Since  the  behavior  of  the  integrand  for  the  chordwise  loading  is  apt  to 
exhibit  similar  characteristics  near  the  surface  edges,  it  is  inferred  that  this 
set  should  also  yield  the  best  approximation  in  three-dimensional,  unsteady, 
compressible  flow.  Note  that  the  number  of  collocation  points  is  not  required 
to  be  the  same  as  the  number  of  integration  points.  As  will  be  seen  later,  it 
is  only  necessary  that  the  total  number  of  downwash  collocation  points  be 
equal  to  or  greater  than  the  number  of  pressure  coefficients  anm»  When  N 
chordwise  integration  stations  are  used,  the  quadrature  used  to  evaluate  the 
inner  integral  of  Equation  51  is  exact  for  integrands  represented  by  a 
polynomial  of  degree  <  2  N  -  1 . 


Hsu  shows  that  an  optimum  set  of  interdigitated  spanwise  collocation 
stations  and  integration  stations  exists  for  evaluation  of  the  outer  integral 
in  Equation  51.  By  reasoning  similar  to  that  used  to  establish  the  chordwise 
collocation  stations,  it  was  established  that  the  optimum  spanwise  collocation 
stations  are 


-r  =  '  COS  MTT  77  ’  r  =  1.  •  •  •  M.  (57) 

It  was  observed  that  the  quadrature  for  the  integral  of  difference  between  the 
actual  and  polynomial  approximation  of  the  spanwise  loading  is  zero  when  the 
actual  loading  is  precisely  represented  by  a  polynomial  of  degree  <  2M  -  1, 
and  the  polynomial  approximation  is  of  degree  =  N  -  1. 
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Then,  by  substitution  of  Equation  55  into  Equation  51, 


/  yd-2 


D  =  —  J  - 5-  G  (x.,  s  ,  cr)  d  <r 

nm  8ir  /  ,  .2  nm  j*  -r*  — 
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where 


2it 
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\07i 


Gnm  ■  WtT.S  am(‘-5k)Un(It)  (ir-<r)2K(xrek,  ir.  2)H  „ 


k  =  1 


Hsu  established  the  form  of  the  Gaussian  quadrature  and  the  spanwise 
integration  stations.  The  difficulties  of  the  singularity' of  the  kernel  at 
?  ~  Sjc  and  the  difficulty  of  differentiation  with  respect  to  or.  (he  uses  the 
steady-state  lifting  line  formula  to  derive  the  form  of  the  quadrature)  are 
avoided  by  removal  of  the  singularity  at  0;  =  sr  and  then  by  an  integration  by 
parts. 


Y/e  first  integrate  by  parts  to  get 


1  JL 

Dl  =  JL  /  -Ail 

nm  8tt  ^ 

which  corresponds  to  Equation  58  in  Reference  9.  The  Gaussian  quadrature 
formula  is  developed  and  shows  that  when  the  number  of  integration  Stations 
is  one  greater  than  the  number  of  collocation  stations,  and  if  they  are  inter- 
digitated  in  the  prescribed  way 
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or2  G  (x  s  -  cr,  co,  M) 
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D*  =  f- 
nm 
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where 
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r  =  1,  2,  .  .  .  M 

Evaluation  of  the  second  term  in  the  brackets  requires  the  observation 
that  the  multiplier  of  in  Equation  48  goes  to  zero  whenever  the  collocation 
point  is  in  the  plane  of  the  doublet  sheet  located  at  the  integration  point. 
Therefore,  the  finite  part  of  the  integral  of  the  K2  term  is  zero,  and  the  entire 
contribution  comes  from  the  Kj  term.  In  this  case,  y(s)  =  y(v)  and 
ri2  =  (sr  -  £)2  =  0. 


It  can  be  shown  that 


K  (x  -  £  .  s  ,  s  ) 
1  *  *  -r’  -r' 


-2,  x  >£ 
0,  x  <  $  . 


Thus,  the  chordwise  integral  which  defines  Gnm  (xj,  JLr>  8r)  *8 
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nm 
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If  the  range  of  integration  is  extended  to  f  =  1  by  making  the  integrand 
zero  for  ?■  >xj,  the  integral  cannot  be  well  approximated  by  a  polynomial 
because  it  has  a  jump  discontinuity  at  |  =  x  j.  To  overcome  this  difficulty, 
we  write 
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G  (x.,  s  ,  s  )  = 
nm  j  — r  —  r 
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The  second  term  here  depends  on  the  integrals 


lj  dl 

(6la) 


x . 

V5j>  =  J 

which  may  be  evaluated  exactly.  W e  have 


This  list  may  be  extended  as  far  as  it  is  needed. 

_The  first  integral  in  Equation  6la  is  considered  as  an  integral  over 
-1  <  1  <1  with  zero  integrand  when  f  >  xj,  and  Equation  55  is  applied.  The 
result  is 
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Equations  59  and  61  are  used  to  evaluate  the  Dj^’s  given  by 
Equation  60.  Equation  50  is  then  used  to  determine  the  pressure  coefficients 
anm  t0  match  the  downwashes  (i.e.  the  w^/U’s)  at  the  collocation  points 
(xj,  sr).  Once  the  pressure  coefficients  are  determined,  the  generalized 
forces  are  computed.  A  polynomial  expression  for  the  i^  modal  deflections 
normal  to  the  surface. 


N  M 


n(i)  (£,<r)  = 


I  I  eV 


v  =  0  u  =  0 

and  Equation  49  are  substituted  into  the  equation 


vu 


S  x 

TIP  TE 


=  /  /  Ap  ^  (i,a_)  d£dcr 


to  get 


Q.. 

ij 


where 


-S  X 

TIP  XLE 


N  M  N  M 

s2TIp(l/2pU2)  lilt  h.'1' 


n  =  0  m  =  0  v  =  0  u  =  0  nm  m  AQnm^ 


(62) 


1  1 


AQ«n ,*‘f  J^4Fn(l)dfd2 


-1  -1 


The  quadrature  formula  given  by  Equation  55  with  L  set  equal  to  N  may  be 
used  to  evaluate  the  inner  integral.  For  evaluation  of  the  outer  integral,  a 
quadrature  formula  with  the  weight  function  y/T~^V},  analogous  to  Equation 
60  for  a  nonsingular  integral,  may  be  used.  This  is  not  practical  for  a  wing 


■j 


29 


with  folded  tip,  for  which  a  different  representation  should  be  used  on  each 
plane  part  of  the  surface.  Here,  the  spanwise  integral  over  one  part  of  the 
wing  may  have  a  square  root  factor  at  one  end  of  the  interval  of  integration, 
or  at  neither  end.  Such  integrals  may  best  be  evaluated  by  a  suitable 
application  of  ordinary  Gaussian  quadrature  with  a  weight  function  of  1.0. 
The  points  and  weights  for  this  quadrature  method  are  not  given  by  simple 
formulas  such  as  Equation  55.  They  are  listed  in  many  places,  (e.  g. 
Reference  5). 
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VI.  APPLICATION  TO  A  WING  WITH  A  FOLDED  TIP 


The  planform  may  be  any  continuous  surface.  However,  the  computer 
program  was  developed  to  treat  either  a  planar  or  nonplanar  planform  of  the 
type  show;-  in  Figure  2.  (Only  one-half  the  planform  is  shown. )  To  facilitate 
modifications  of  the  program,  comment  cards  are  placed  throughout  the 
program  to  indicate  where  changes  may  be  made  to  handle  other  nonplanar 
surfaces  like  that  of  the  Paraglider. 

In  the  application  of  the  kernel  function  method  to  the  planform  shown 
in  Figure  2,  the  computer  program  calculates  from  the  equations  of  the 
leading  and  trailing  edges  the  collocation  and  integration  points  for  which 
the  integrands  of  Equation  51  must  be  evaluated.  For  demonstration 
purposes,  we  calculate  the  collocation  and  integration  points  for  values  of 
L  =  4,  N  =  6,  and  M  =  10  and  show  the  results  in  Figure  3. 

*LE  =  s  tanX  LE(  *LETIp  =  .  tan  X  LE  +  (s  -  a  )“">•  LET1p 

£j  JL 

*TE  *  2V  3  TE.  XTETIP  =  2V  SF ,*an  ATE+  *S  '  SF,  )taaK  TETIF 

Jb 


b(s) 

=  1/2  *xte  "  x: 

X 

m 

=  l/2(xLE  +  x, 

X 

=  x  +  b(s)  x 
m 

e 

=  em  +  b(o-)? 

cos  2 jyj  j  j  =  1 , 2, . .  .  N. 


2k  -  1 


+  j  v,  k  =  1 , 2, . .  .  L. 
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Figure  2.  Planar  Wing  With  Planar  Symmetrically  Folded  Tips 


•  INTEGRATION  STATIONS  (L- 4) 

•  CHORDWISE  COLLATION  STATIONS  (N  -  6) 
AND  SPANWISE  INTEGRATION  STATIONS 

+  DOWNWASH  COLLOCATION  POINTS  (M  -  10) 


Figure  3.  An  Optimum  Set  of  Collocation  and  Integration  Points 
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rrnr  ,  _ 

-  cos  ~ f  i  m  =  1,2,..  .„M. 


cos  2(iyt  4.  ir>  P  =  1 »  2, .  .  .4  (M  +  1). 


We  have  also  constructed  a  table  of  equations  (Table  1),  for  the 
functions  yQ,  zQ,  r^  Tj,  and  T2  for  use  in  the  expression  for  the  kernel, 
Equations  48  through  48c.  In  Table  1,  the  subscript  Fj^  indicates  the  fold 
line. 


A  difficulty  is  encountered  in  the  spanwi.se  integration  because  the 
kernel  function  has  a  finite  discontinuity  at  the  fold  line.  For  example,  note 
in  Table  1  the  change  in  Tj  and  T2  for  receiving  or  collocation  points  on  the 
wing  as  the  sending  or  integration  points  shift  from  the  port  tip  to  the  wing 
and  from  the  wing  to  the  starboard  tip. 

If  we  consider  the  kernel  as  a  function  of  £  and  £  for  fixed  values  of  x  and  a, 

q{£  ,  _a)  =  K  (x  -  £  ,  s,£), 

this  function  may  be  broken  up  into  a  simple  discontinuous  part  g**  and  a 
part  g*  which  is  continuous  across  the  fold  lines: 


To  do  this,  define 


g(S.£)  =  g*(|  »£.)  +  g**(£  •  £) 


g(6  »£F  +)  "  g(€  .£F  -).  £  >  £F 
L  L  L 


g  **(?,£)  =  0, 


-£F  <*F 
JL  .L. 


g($.-£F  -)  -  g (£.-£F  +).  <L  <  -SL~ 

L  L  L 


and  then  define  g*  by  Equation  63. 


More  explicitly,  for  0-  >  o-  ,  define 

~FL 


£f  =  b(£  )  £  +T  )  +  ZwiZr  ) 


2  b  LE  F.  '  '  TE^F  ' 
*  A-  JL 
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are  given  by  formulas 


f  J.  F  m  , 

U  =  /  »  1  -  (T  0[  a£ 


1  IT  .  -1 

u„  =  —  —  -  sin  £_  +  jjt  t' 

0  2  2  fl  fl 


1  -  J£. 


1  ,,  2  3/2 

“l  3  11  ‘  “F^ 


In  terms  of  these  constants 


mn  8ir 


A 

2  m  f 

-  jr  a;  J 


fn(S)g*(S,<r)d£dcr 


+  isrv  / fn<«>  ksVe)  +  M)mg“(?)]  o? 


The  last  integral  in  this  formula  is  evaluated  by  Equation  55. 

In  the  evaluation  of  Ojj,  given  by  Equation  62  for  the  plane  case,  it  is 
assumed  that  modes  numbers  i  and  j  are  either  both  symmetric  in  it,  or  both 
antisymmetric.  Then  the  contribution  to  the  integral  for  cr  <  0  is  the  same 
as  the  contribution  forjr  >  0.  Let  the  deflection  in  the  i^  mode  be  given  by 


N  M 


2  2  c„'  0  <  ff  <  21 F 


n(i)(e,2l)  = 


v  =  0  p  =  0 


N  M 


I  2 


v  =  0  u  =  0 
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Then 
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ij 


in  which 
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Note  that  the  inner  integrals  are  the  integrals  of  polynomials  in  £  multiplied 
by  v(l  -  £)/(l  +  £),  by  virtue  of  the  relation 

£  =  b(ff)g  +  j  [$le(jl)  +  tTE(s;)] 


Hence,  the  inner  integral  may  be  evaluated  exactly  by  either  Equation  52a  or 
Equation  55  if  enough  points  are  used  in  the  formulas.  For  the  limits  v  <  5, 
n  S  4  used  in  the  computer  program,  six  points  are  sufficient  for  Equation 
52a,  five  points  for  Equation  55.  Equation  52a  was  used  with  K  =  6.  (.This 
choice  was  made  arbitrarily;  it  would  be  just  as  good  to  use  Equation  55.  ) 

In  the  integrations  over  jr  ,  six-point  Gaussian  integration  with  weight 
function  1.  0  was  used.  The  basic  formula  is 


1  6 

/  f(v)  dv  =  X  (66) 

0  £=  1 

exact  for  f(v)  a  polynomial  of  degree  at  most  11.  The  constants  occurring 
in  this  formula  are  given  in  the  subroutine  FORCE.  -They  may  be  derived 
from  a  table  given  by  Scarborough  (Reference  6,  p.  148). 


In  ,  Equation  66  is  applied  by  putting 

nmvn 


<r  =  £  v 

fl 


The  resulting  expression  is 
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nmvp. 
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V  ,  r  ~Z  m  +  to. 

Z  hi  x/l  -z*  £> 


-i  1 
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in  which 


“I 


^k(<rt)  b^£^k+  2  +  ^TE^l 


(2) 

In  I  .  the  transformation 

nm  vu. 


£  =  l-  (l-£_)v 


was  used.  This  makes  the  v-integrand  behave  like  a  polynomial  at  the  ends 
of  the  interval.  The  resulting  formula  is 


t(2)  _  / 1  1 

nmv(i  -  3  •-Fl)  =  j 
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DESCRIPTION  OF  THE  COMPUTER  PROGRAM 


A  functional  diagram  of  the  computer  program  is  given  in  Figure  4. 
With  the  exception  of  two  subroutines  named  MSIMEC  and  MSIMER,  all  of 
the  programs  are  written  in  Fortran  IV.  These  two  subroutines,  written  in 
machine  language,  are  used  for  complex  and  real  matrix  inversion, 
respectively.  There  are  certain  limitations  related  to  the  various  other 
subprograms,  which  are  listed  below. 


Subprogram 


Limitations 


Subroutine  Data 

NCC 

The  number  of  chordwise  collocation  stations  must 

be  £  10. 

NCS 

The  number  of  spanwise  collocation  stations  must 

be  £  9. 

NDATA 

The  number  of  sets  of  data  must  be  £  10. 

N 

The  number  of  chordwise  pressure  modes  £  5. 

M 

The  number  of  spanwise  pressure  modes  <  5. 

Subroutine  Zen 


MODES 


NPTS 


This  is  the  number  of  modes  used  in  the  calculation 
of  generalized  forces,  and  must  be  £  10. 

For  a  planar  wing,  this  is  the  number  of  points  at 
which  the  deflection  is  given  in  the  horizontal 
surface  and  must  be  £  66.  For  a  nonplanar  wing, 
there  must  be  £  66  points  for  the  deflections  in  the 
horizontal  surface  and  £  66  points  for  the  deflections 
in  the  vertical  surface. 


*2  O 


The  first  set  of  data  is  read.  This  includes  Items  1-14  In  the 
data  array,  the  Mach  numbers  and  frequencies,  item  41 ,  and 
perhaps  some  of  the  mode  data  that  will  be  needed  later. 


DATA1 


Subroutine  DATA1  performs  most  of  the  calculations  that  ore 
independent  of  Mach  number  and  frequency,  such  as  setting  up 
the  arrays  of  collocation  and  integration  stations. 


This  is  the  beginning  of  a  loop  that  is  gone  through  for  each  com¬ 
bination  of  Mach  number  and  frequency  specified.  Subroutine 
DATA2  sets  up  the  matrix  (D'^)  of  Equation  50. 


In  subroutine  ZEN,  any  necessary  data  for  the  oscillation  modes 
is  read  in.  The  left  side  of  Equctlon  50  Is  set  up  for  eoch  mode. 


<3H 


MSIMEC 


Equation  50  Is  solved  for  the  pressure  coefficients  for  each  mode 
by  MSIMEC. 


F0RCE 


Subroutine  F0RCE  computes  the  generalized  force  matrix.  The 
program  returns  to  step  3  until  the  Mach  numbers  and  frequencies 
are  exhausted,  then  to  step  1 . 


(*)  RETURN  TO  BEGINNING  FOR  NEXT  CASE 

(«*)  LOOP  ON  NUMBER  OF  SETS  OF  MACH  NUMBERS 
AND  FREQUENCIES 


Figure  4.  Functional  Flow  Diagram — Main  Program 


USE  OF  THE  COMPUTER  PROGRAM 

The  following  rules  apply  to  optimum  use  of  the  computer  program': 

1.  Before  attempting  to  use  the  computer  program  compute  the 
coefficients  of  polynomials  of  minimum  order  in  x  and  y  that 
will  adequately  represent  the  modal  deflection  distributions. 
Least-square  fitted  surfaces  are  very  useful  for  this  purpose, 
since  weighting  factors  may  be  used  to  obtain  a  better  fit  in 
special  regions. 

2.  Set  the  number  of  chordwise  collocation  points  M  equal  to  one 
plus  the  highest  of  the  orders  in  x;  and,  unless  there  is  special 
reason  to  reduce  the  number  of  chojrdwise  integration  points,  set 
L  equal  to  M. 

3.  Set  the  number  of  spanwise  collocation  points  R  equal  to  one 
plus  the  highest  of  the  orders  in  y.  This  establishes  the  number 
of  pressure  coefficients  anm  at  M  x  R.  Their  values  are 
computed  by  matching  exactly  the  downwashes  at  file  M  x  R  span- 
wise  collocation  points. 

The  Input  .data  is  read  by  the  subroutine  DATRD.  Use  of  this  sub¬ 
routine  requires  that,  on  each  data  card,  the  first  72  columns  are  six 
fields  of  width  12,  as  indicated  on  the  sample  date  sheets  (Figure  5). 

The  first  field  contains  an  integer  giving  the  location  in  the  data  array  in 
which  the  number  in  the  second  field  is  to  be  scored.  The  numbers  in  the 
remaining  fields  are  stored  in  consecutive  locations.  If  a  field  is  blank,  the 
corresponding  location  in  the  data  array  is  unchanged.  DATRD  reads  any 
number  of  cards.  A  minus  sign  in  column  1  indicates  the  last  card  to  be 
read;  if  this  minus  sign  is  not  present,  DATRD  continues  with  the  next  card. 
The  storage  locations  of  the  data  on  a  card  are  not  affected  by  the  sign  in 
column  1.  All  floating  point  numbers  must  be  written  with  decimal  points. 
All  integers  must  be  at  the  right  of  their  fields. 

The  data  array  is  set  up  as  follows: 


1. 

N 

The  number  of  chordwite  pressure  motfes 

2. 

M 

The  number  of  spanwise  pressure  modes 

3. 

NCC 

The  number  of  chordwise  collocation  points 

4. 

NCS 

The  number  of  spanwise  collocation  points 

6. 

NDATA 

The  number  of  sets  of  values  of  Mach  number  and 
frequency  to  be  used. 

6. 

NSYM 

Indicator  for  symmetric  (NSYM  *  +1)  or 
antisymmetric  (NSYM  =  -1)  modes  of  oscillation. 
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7. 

SFOLD 

Distance  spanwise  from  wing  center  line  to  fold  line 

8. 

S'flP 

Semispan 

9. 

BO 

One-half  of  the  root  chord 

10. 

ALFA1 

The  fold  angle  (in  degrees) 

t  U. 

i 

i 

^LE 

The  sweep  angle  of  the  leading  edge  (in  degrees) 

|  12. 

t 

\TE 

13. 

^LEup 

14‘ 

^tetip 

21-30. 

Values  of  Mach  number. 

NDATA  of  these  must  be  entered. 

31-40, 

Values  of  reduced  frequency 
wBOAJ. 

NDATA  of  these  must  be  entered. 

41. 

NMOD 

The  number  of  modes 

42. 

JD 

Indicator  for  the  type  of  input  data  for  a  mode.  If 
JD=1,  the  deflections  are  given  at  a  set  of  points  on 
the  wing  (and  tip).  If  JD=2,  the  coefficients  of 
polynomials  for  the  deflection  of  the  wing  and  tip 
are  given.  If  JD=0,  the  current  mode  and  subsequen 
modes  are  not  given  by  data.  They  are  the  same  as 
the  corresponding  modes  which  were  used  for  the 
previous  frequency  and  Mach  number. 

43. 

NPTSW 

Number  of  points  on  the  wing  at  which  deflections 
are  given  (used  if  JD=1). 

44. 

NPTST 

Number  of  points  on  tip  at  which  deflections  are 
given. 

51-71. 

Deflection  coefficients  on  the  wing. 

The  coefficients  are  stored  as  follows: 

51  52  53  54  55  56  57  58  59 

a.„  +  a,„  x  +  a_  x2  +  a.,  x3  +  x4  +  a,,  x5  +  y(a.,+a,,x+a„  x2 
00  10  20  30  40  50  01  11  21 


60  61 


62  63  64  65 


+  a„,  x3+  a„  x4)  +  y2(a-„+a,„x+a  x2+a„,x3) 

31  41  '  7  ^  02  12  22  32  ' 


66  67  68 


69  70  71 


+  y3(a„„+a,„x+a„„x2)  +  y4(a  +a,  x)  +  a„.y5 
1  '  03  13  23  '  1  v  04  14  '  057 
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76-96.  Deflection  coefficients  on  the  tip.  Same  storage  rule  as  above,  except  all  locations 

are  increased  by  25. 

98.  Indicator  that  no  more  modes  are 
to  be  read  after  the  present  one. 

(For  current  and  subsequent 
frequencies  and  Mach  numbers.) 

101-299.  Deflection  data  at  points  on  the 
wing,  in  the  order  x^  s  ,  n  , 

WetC* 

301-499.  Deflection  data  at  points  on  the  tip. 

Items  1-41  must  be  read  in  the  first  set  of  data.  There  may  be  an 
additional  set  of  data  for  each  mode,  for  each  Mach  number  and  frequency 
case.  After  the  indicator  DAT(98)  has  been  given  a  non-zero  value,  no  more 
deflection  data  will  be  read.  After  JD  has  been  given  the  value  zero,  no 
data  will  be  read  for  the  higher  numbered  modes. 

The  data  in  Figure  5  is  for  a  60°  triangular  wing  folded  at  75  percent 
semispan,  at  an  angle  of  30°.  The  root  chord  is  5.  0  feet,  making  BO  =  2.  5. 
Three  modes:  plunge,  pitch,  and  a  third  nonrigid  mode  are  considered.  The 
Mach  number  is  0.7,  and  six  frequencies :  10,  20,  30,  40,  50,  and  60  cps 
are  used.  The  speed  of  sound  is  taken  to  be  1000  ft.  /sec.  This  gives  reduced 
frequencies  of  0.157,  0.314,  0.471,  0.628,  0.785,  and  0.942. 

Three  spanwise  and  three  chordwise  pressure  modes,  six  chordwise 
and  eight  spanwise  collocation  stations  are  specified. 

In  the  set  of  cards  numbered  15  -  30,  which  give  deflection  data,  some 
of  the  cards  have  been  omitted.  Otherwise,  this  is  a  complete  set  of  data  for 
a  computer  run. 
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FORTRAN  FIXED  10  DIGIT  DECIMAL  DATA 

DECK  NO. _  PROGRAMMER  _  OATE_  PAOE_JL_— of.  5  JOB  NO. 

NUMBER  IDENTIFICATION  DESCRIPTION  DO  NOT  KEY  PUNCH 
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Coefficients  of  the  deflection  polynomial  on  the  wing,  which  has  the 
constant  value  1.0  In  this  mode. 


indicates  last  card  for  first  mode 


FORTRAN  FIXED  10  DIGIT  DECIMAL  DATA 

DECK  NO...  . .  .  PROGRAMMER  -  DATE  _  PAGE  4  of  5  JOB  NO. 

NUMBER  IDENTIFICATION  DESCRIPTION  DO  NOT  KEY  PUNCH 


Figure  5.  Sample  Data  Sheets  (Sheet  4  of  5) 


VII.  RESULTS 


The  computer  program  was  applied  to  the  rigid  modes  of  two  wings. 

ASPECT  RATIO  2.  0  RECTANGULAR  WING  FOLDED  AT  80  PERCENT 
SEMISPAN 

CLq4s  plotted  as  a  function  of  Mach  number  in  Figure  6.  The 
dashed  curve  is  a  plot  of  the  approximation 


CL„  = 


2ttA.  R. 


2  4  +  (A.  R.  )2  (1-M2) 


where  A.R.  is  the  aspect  ratio.  This  is  formula  (6-31)  of  Reference  5,  for 
the  case  of  a  rectangular  wing. 

TRIANGULAR  WING  WITH  FOLDED  TIPS 

The  configuration  used  was  a  triangular  wing  with  a  sweep  angle  of 
65  degrees,  folded  at  60  percent  semispan. 

Figures  7  and  8  show  CL^and  Cmaas  functions  of  fold  angle 
(at  M  =  0.  8),  and  as  functions  of  Mach  number.  Figure  9  is  a  plot  of 
unsteady  generalized  forces  for  rigid  oscillations  of  the  wing  in  the  pitching 
mode. 
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Vill.  CONCLUSIONS  AND  RECOMMENDATIONS 


The  results  obtained  by  the  nonplanar  kernal  function  method  show 
the  expected  trend  with  increasing  fold  angle,  which  agrees  with  the 
observed  evperimental  trend.  For  zero  fold  angle,  the  method  reduces  to 
that  already  used  by  Hsu  (Reference  4). 

Possible  extensions  of  the  method  include  the  treatment  of  more 
general  configurations,  or  of  other  specific  configurations,  such  as  the  T- 
tail.  Also,  any  generalizations  proposed  for  the  planar  case  should  be 
considered  here,  such  as  the  problem  of  a  nonplanar  wing  with  a  control 
surface. 

The  formula  that  was  used  for  the  kernel  function  (page  19)  should  be 
useful  in  any  future  developments  using  the  three-dimensional  kernel  function. 
The  previously  available  formula  was  much  longer.  A  special  case  of  that 
formula  is  given  in  Reference  12.  The  use  of  the  simplified  kernel  function, 
together  with  Hsu's  method  of  integration,  results  in  greatly  reduced  com¬ 
puter  running  times.  This  makes  it  practical  to  use  the  kernel  function 
method  as  a  tool  in  the  preliminary  analysis  of  new  wing  configurations. 
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